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Multi-channel spin-orbit quantum wires, when subjected to a magnetic field and proximity cou- 
pled to s-wave superconductor, may support Majorana states. We study what happens to these 
systems in the presence of disorder. Inspired by the widely established theoretical methods of 
mesoscopic superconductivity, we develop a la Eilenberger a quasiclassical approach to topological 
nanowires valid in the limit of strong spin-orbit coupling. We find that the "Majorana number" M, 
distinguishing between the state with Majorana fermion (symmetry class B) and no Majorana (class 
D), is given by the product of two Pfaffians of gapped quasiclassical Green's functions fixed by right 
and left terminals connected to the wire. A numerical solution of the Eilenberger equations reveals 
that the class D disordered quantum wires are prone to the formation of the zero-energy anomaly 
(class D impurity spectral peak) in the local density of states which shares the key features of the 
Majorana peak. In this way we confirm the robustness of our previous conclusions [Phys. Rev. Lett. 
109, 227005 (2012)] on a more restrictive system setup. Generally speaking, we find that the qua- 
siclassical approach provides a highly efficient means to address disordered class D superconductors 
both in the presence and absence of topological structures. 

PACS numbers: 74.78.Na, 71.23.-k, 73.63.Nm, 74.45. +c 



I. INTRODUCTION 



Semiconductor quantum wires proximity coupled to 
a conventional superconductor and subject to a mag- 
netic field may support Majorana fermion edge states^. 
Building on relatively conventional device technology, 
this proposed realization of the otherwise evasive Ma- 
jorana fermion has triggered a wave of theoretical and 
experimental activity, which culminated in the recent re- 
port of a successful observation by several experimental 
groups®^. In these experiments, evidence for the pres- 
ence of a Majorana is drawn from the observation of a 
zero bias peak in the tunneling conductance into the wire. 
While the observed signal appears to be naturally ex- 
plained in terms of a Majorana end state, two of us have 
pointed out that a midgap peak might be generated by 
an unrelated mechanism^ in the presence of even very 
moderate amounts of disorder, the semiconductor wire 
supports a zero energy "spectral peak" (an accumulation 
of spectral weight at zero energy) which resembles the 
Majorana peak in practically all relevant aspects. Specif- 
ically, it is (i) rigidly locked to zero energy, (ii) is of nar- 
row width of 0(6), where 8 is the single particle level 
spacing, (hi) carries integrated spectral weight 0(1), and 
(iv) relies on parametric conditions (with regard to spin 
orbit interaction, proximity coupling, magnetic field, and 
chemical potential) identical to those required by Majo- 
rana state formation. What makes the spectral peak dis- 
tinct from the Majorana is that it relies on the presence 
of a moderate amount of disorder, viz. impurity scatter- 
ing strong enough to couple neighboring single particle 
Andreev levels. Besides, the spectral peak is vulnerable 
to temperature induced dephasing. While this marks a 
difference to the robust Majorana state, the reported ex- 
perimental data does show strong sensitivity to temper- 
ature, which may either be due to an intrinsic sensitivity 



of the peak, or due to a temperature induced diminishing 
of the measurement sensitivity, or both. In either case, 
the situation looks inconclusive in this regard. 

Generally speaking, the results of Ref. 6, as well as 
those of Refs. [7| and |8] suggest that the observation of a 
midgap anomaly in the tunneling conductance might be 
due to either mechanism, disorder peak, Majorana peak, 
or a superposition of the two, and this calls for further 
research. 

Our previous study was based on an analytically 
tractable idealization of a semiconductor quantum wire 
subject to a magnetic field sweep. In the present paper, 
we will explore the role of disorder scattering within a 
model much closer to the experimental setup. The price 
to be payed for this more realistic description is that a 
fully analytic treatment is out of the question. Instead, 
we will employ a semi-analytic approach based on the 
formalism of qu asicla ssical Green functions. Introduced 
in the late sixties^^, the latter has become an indispens- 
able tool in mesoscopic superconductivity^^ and quan- 
tum transport in general. We here argue that quasiclas- 
sical methods are, in fact, tailor made to the modeling of 
Majorana quantum wires (or, more generally, quasi one- 
dimensional topological superconductors.) Specifically, 
we will show that in the problem at hand, the quasiclas- 
sical "approximation" is actually very mild. Further, the 
quasiclassical Green's function affords a convenient de- 
scription of the topological signatures of the system in 
terms of Pfaffians. Finally, the approach can be applied 
to systems for a given realization of the disorder, and at 
numerical cost much lower than that of exact diagonal- 
ization approaches. As a result, we will be in a position 
to accurately describe local spectral properties within a 
reasonably realistic model of a topological multi-channel 
superconductor. As we are going to discuss below, our 
findings support the principal statement made in Ref. 6 
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The rest of the paper is organized as follows. In sec- 
tion [TlJ we discuss the principal role played by disorder 
in the system, the idea of the quasiclassical approach, 
and its main results. In section |TH| we specify our model 
system, the quasi-classical approach is introduced in sec- 
tion |IV[ and in [V] we discuss the numerical solution of 
the quasiclassical equations. We conclude in section |VT] 
A number of technical details are relegated to Appen- 
dices. 



II. QUALITATIVE DISCUSSION AND RESULTS 

A schematic of the device currently under experimen- 
tal investigation is shown in Fig. [T] A semiconductor 
quantum wire subjected to strong spin orbit interaction 
is brought in contact to a superconductor (S), and, via a 
tunnel barrier (T) to a normal metal lead (N). The appli- 
cation of a small excess voltage, V, to the latter induces 
a tunnel current into the central region. The differen- 
tial conductance dl/dV probes the (tunneling) density 
of states at an energy V (units e = h = c = 1 through- 
out) relative to the systems chemical potential, /i. The 
physics we are interested in is contained in a band center 
(V = 0) anomaly in that quantity. 

In a manner to be discussed in more detail below, one 
contribution to the band center density of states is pro- 
vided by a Majorana bound state localized at the tunnel 
barrier. The second, spectral contribution is generated 
by a conspiracy all other low lying quasiparticle states in 
the system. Technically, these are Andreev bound states 
forming at energies ±Ej in the region between the tun- 
nel barrier and the superconductor. The number of these 
states increases with the extension of the wire - a few 
hundred nanometers in the experiment - and the number 
of transverse channels below the chemical potential. In 
the presence of disorder, the entity of these states defines 
an effective "quantum dot". The proximity of a super- 
conductor, and the breaking of both spin rotation and 
time reversal symmetry imply that the system belongs 
to symmetry class 

The symmetry of class D random systems implies a 
clustering of levels at zero energy. Loosely speaking, the 
conventional level repulsion of random spectra turns into 
a zero energy level attraction. On a resolution limited 
to scales of order of the mean level spacing, the zero 
energy density of states (DoS) is enhanced by a factor 
of two relative to the mean background, i.e. it shows 
a peak. This phenomenon manifests itself at relatively 
small sample-to-sample fluctuations, i.e. the peak is a 
sample specific effect. In passing we note that the weak 
anti-localization phenomenon discussed in Ref. [5] rests 
on the same principal mechanism of midgap quantum 
interference. 

Below, we will explore the phenomenon of spectral 
peak formation in a setting modelled to closely mimic the 
"experimental reality" . To be more specific, we consider 
a semiconductor wire supporting a number of TV > 1 





















TV 




T 































T L 

— ^ 

i A 



I l V 

FIG. 1. Disordered spin-orbit nanowire, subjected to colinear 
magnetic field, is proximity coupled to the s-wave supercon- 
ductor (S) and terminated by the tunneling barrier (T) at one 
of its ends. A sketch below shows the profile of the induced 
superconducting gap A(x) and gate induced potential V(x) 
defining the tunnel barrier. Andreev bound states (ABS) are 
depicted by dotted lines. 



transverse channels below the chemical potential. We 
assume a value of the chemical potential such that the 
highest lying of these channels is "topological" (chemical 
potential falling into the gap opened by the simultaneous 
presence of order parameter and magnetic field.) 

To quantitatively describe this phenomenon, we gen- 
eralize the quasiclassical Green function approach to the 
present setting. Indeed, the phenomena we are interested 
in manifest themselves on length scales large in compar- 
ison to the Fermi wavelength, yet smaller than the rele- 
vant coherence length, and this makes them suitable to 
quasiclassical treatment. The quasiclassical theory of the 
present paper is formulated in terms of the Eilenberger 
function Q(x]e), which is a position and energy depen- 
dent matrix of size 8TV x 8TV. The factor of 8 = 2x2x2 
accounts for spin, chiral (left /right modes) and particle- 
hole degrees of freedom. It will turn out that the struc- 
ture of the theory is most transparently exposed in the 
so-called Majorana basis, where the Eilenberger function 
becomes a real antisymmetric matrix. The Pfaffian of 
that matrix, Pf(Q) will be seen to assume two values 
(±1), which locally (in space) identify the I12 invariant 
of the underlying one-dimensional class D superconduc- 
tor, in dependence on system parameters: for values of 
the magnetic field smaller or larger than a critical field, 
B c = \J A 2 + /1 2 , where [i and A are chemical poten- 
tial and bulk order parameter, the invariant is trivial or 
non-trivial, resp. In the latter case, the wire supports a 
Majorana state at the barrier, in the former it doesn't. 

Fig. [2] shows a profile of the local DoS (LDoS) com- 
puted at the left end of the wire for typical system pa- 
rameters as detailed below. The green dashed and red 
solid curve correspond to a situation without and with 
Majorana state. Within our model, the states in the wire 
carry a narrow width, ~ gr^, reflecting the possibility of 
decay through the tunnel barrier into the adjacent lead. 
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FIG. 2. Disorder averaged local density of states (LDoS) at 
the left end of the spin-orbit quantum wire sketched in Fig. [I] 
(in units v = l/2nv). The number of occupied bands N = 2 
corresponds to four transport channels. Parameters are: (i) 
red solid line (A4 = —1), B — 2.66A, jj, = A — topological 
phase; (ii) green dashed line (M = +1), B = 0.25A, (j, = 0.5A 
— trivial phase. The wire length L — 4v/A, dimensionless 
strength of disorder 7„/v = 0.16A, which translates into the 
mean free path I = 0.4L. Tunneling rate T = 5 ■ 10~ 2 A. 
Velocities in two bands were taken to be equal, Vi = V2 = v. 
The inset shows profiles of the DoS resulting from random 
matrix theory. 



(Here, 8 is the mean level spacing of the wire, and gr < 1 
the barrier tunneling conductance.) This broadening ac- 
counts for the finite width of the Majorana peak in the 
topological regime. Loosely speaking, the negative shoul- 
ders observable next to the center peak are due to the 
repulsion of adjacent levels of the center Majorana level. 
A more substantial explanation is as follows: for an odd 
parity of the total number of levels — a signature of the 
topological regimes — the disordered quantum system 
falls into symmetry class B, rather than D. (Class B is 
the designation for a system with the same symmetries 
as a class D system, yet odd number of levels.) A uni- 
versal signature of class B is a negative spectral peak at 
zero energy (the negative of the positive class D peak), 
superimposed on a single (5-peak (the Majorana). The 
joint signature of these two structures is seen in the solid 
curve in Fig. [2] 

At resolutions limited to values ~ S, the superposition 
of the Majorana and the class B peak looks next to in- 
distinguishable from the class D peak (dashed), and this 
similarity of unrelated structures might interfere with the 
unambiguous observation of the Majorana by tunneling 
spectroscopy. Indeed, the differential tunneling conduc- 
tance monitored in experiment, 
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FIG. 3. The average LDoS (solid green line) and the square 
root of its second (reducible) moment (^£(e)} 1,/2 (dashed blue 
line) at the left end of the spin-orbit quantum wire in the 
trivial phase (A4 = +1). System parameters are listed in 
Fig. [2] The inset shows profiles of the mean DoS and the 
square root of the two level correlation function resulting from 
random matrix theory. 



is essentially^ determined by the local density of states, 
i.e. the structures shown in Fig. [2] are expected to reflect 
directly in the measured signal. (Here, fp is the Fermi 
distribution, vl is the DoS at the left barrier, v is the 
DoS in the single chiral channel per unit length, and V 
the applied voltage.) 

The profiles of the curves shown in Fig. [2] were com- 
puted for a two-channel quantum wire at a mean free 
path I ~ L of the order of the system size. We are ad- 
dressing a system at the interface between the ballistic 
and the localized regime. In view of these system param- 
eters it is remarkable that the DoS profiles in Fig. [2] show 
striking similarity to the average DoS of a class D and 
B random matrix modePA For comparison the average 
DoS of a class B and D random matrix Hamiltonian is 
shown as an inset in Fig. [2] The similarity of the results 
indicates that the system of subgap states in our system 
behaves as if it formed an effective chaotic quantum dot 
localized in the vicinity of the left system boundary (right 
to the tunnel barrier). Taking into account spin, chiral 
and channel quantum numbers, the mean level spacing 
in such dot is given by 6 ~ irv/2NL. In the presence 
of magnetic field the BCS gap in the superconducting re- 
gion of the wire reads e_ = \B — \J A 2 + /i 2 |. The number 
of subgap Andreev levels forming the effective dot is thus 
given by AWls - 2e-/5. 

The profiles shown in Fig. [2] are ensemble averages, 
(vl) of the LDoS, v^, where the sampling was over ~ 500 
randomly chosen impurity configurations. To demon- 
strate the weakness of fluctuations, Fig. [3] compares the 
average LDoS (solid line) to the "typical" LDoS, i.e. the 
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FIG. 4. The sample specific LDoS in the trivial phase without 
Majorana state (M = +1) for two different disorder realiza- 
tions. Tunneling rate T — 0.05 A, other system parameters are 
the same as in Fig. [2] Curves demonstrate two typical sce- 
narios: (i) two conjugate Andreev bound states ±e m in lying 
close to Fermi energy and having energy splitting ~ F (dotted 
green line); (ii) particle and hole states have merged into a 
single zero-energy peak of width F and can not be resolved 
by tunnel spectroscopy (solid magenta line). For the chosen 
set of parameters the mean level spacing 5 — 0.2A and the 
gap in the S region e_ = 0.87A. Thus one has approximately 
iVieveis — 8 random Andreev levels. 



average \f\v\). The relatively minor deviation between 
average and typical DoS demonstrates that the standard 
deviation 



5v{e) = ((v{e)-{v{e)))' 



2\ 1/2 



(2) 



characterizing the strength of mesoscopic fluctuations is 
relatively small. 

The weakness of fluctuations implies that the disor- 
der peak is a realization specific phenomenon. This 
is demonstrated explicitly in Fig. [4j where two un- 
averaged DoS profiles individual disorder configurations 
are shown. The data is for a "non-topological" system, 
no non-degenerate zero energy level is present. Depend- 
ing on whether the pair of lowest lying Andreev states 
(+emini — emin) exceeds the level broadening T, the DoS 
enhancement will assume the form of a single peak (ma- 
genta solid), or a split peak (green dashed). In either 
case, an excess DoS of integrated spectral weight ~ 1 is 
present. 

We finally note that both the Majorana peak and the D 
spectral peak crucially rely on the presence of a magnetic 
field. In the absence of a field time-reversal symmetry is 
restored, and the wire turns into a member of symmetry 
class DHL Such systems display a DoS depletion at zero 
energy, rather than a peak. In other words, the spectral 
peak discussed here will disappear along with the Majo- 



rana resonance. 

In the rest of the paper we discuss how the results sum- 
marized here were obtained by a combination of quasi- 
classical and numerical methods. 



III. MODEL OF THE MAJORANA NANOWIRE 



We consider a multi-band quantum wire of width L z , 
subject to Rashba spin-orbit coupling, proximity cou- 
plingto an s-wave superconductor, and to a magnetic 
fields. We choose coordinates such that the wire lies 
along the x-axis, parallel to the magnetic field, the y-axis 
is perpendicular to the surface of the superconductor, and 
the spin orbit field is pointing along the z-axis. 

In the rest of this section we will specify the 
Bogoliubov-de-Gennes (BdG) Hamiltonian describing 
this system in the presence of disorder (section III A I. 



We will then linearize the electron spectrum in the limit 
of strong spin-orbit coupling thereby introducing a de- 
scription in terms of right and left one-dimensional chiral 
fermions (section IIIB I, and finally transform the Hamil- 
tonian into Majorana basis (section III C I. Experts in the 
description of topological quantum wires may proceed di- 
rectly to the end of the section. 



A. Bogoliubov-de-Gennes Hamiltonian 

Introducing a four component spinor = 
(ip^, ip±, ip-f, Vj,) in the product of spin and particle- 
hole spaces, the BdG Hamiltonian H describing the 
system in the xz-plane reads 



h Q + W 
—is y A 



is y A* 



W 1 



with 
h = 



^(x, z)dx dz, 
(3) 



-{dl + dl)/2m - n{x) + B(x)s x - ia{s z d x - s x d z ). 



Here, A = A(x) is the proximity amplitude induced the 
superconductor, B = ^g/isH, where H is the external 
field, we have taken into account the transverse momen- 
tum (—id z ) in the Rashba term, and W{x) is the random 
disorder Hamiltonian. The Pauli matrices s operate in 
spin space. 

We proceed by introducing a system of transverse wave 
functions, {^(z)}, which leads to a linear decomposition 
of the Grassmann fields (a =f, J.) as 



Defining 



&=5>*(*)*& n) 0r). (4) 



^») = (^,Vf ) ,< ) ! ^ n) ) T , 
the Hamiltonian then takes the form 



(5) 
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where the one-dimensional Hamitonian h acts in the 



-th band as 



;(«) 



<9 2 /2m + _ ^ + B^)^ - ias*^. (7) 



In the case of an ideal waveguide the transverse wavefunc- 
tions are = yj2/L z sm(nirz/L z ), so that fi z n ^ = 

7r 2 n 2 /(2mL 2 ) and the matrix elements of the spin-orbit 
interaction read 



2ia 



(1 - (-l) n+m ) s x 



(8) 

Let us now assume a thin wire, L z < ^ so = h/(ma), 
where l so is the spin-orbit length. In this case the ma- 



trix elements 



•C Hz"' can be treated as perturba- 
tions. To this end we introduce a unitary transforma- 
tion U, which brings the high energy part of the Hamil- 

tonian, (fii '<5„ m + him), to diagonal form. Due to 
the weakness of the perturbation, the transformation is 
close to unity, U = exp(iX) ~ 1 + iX, with generators 
X nm — 0(L z /l so ). To first-order perturbation theory 
their explicit form reads 



X 



■us.o. 



(m) 



fl 



X ,, 



0. 



(9) 



The transformation W — > U^> and \P 
the approximate Hamiltonian 



generates 



H~ L - j¥ n \x) 



w n 



-is y A6 nm -(5A)* nm -(h^yS,. 



(5A) r 

V T - 



w 1 



¥ m \x)dx. 



(10) 



Here we have redefined our notation for the disorder 
potential W — > UWW , and introduced a small correc- 
tion to the quasiparticle Hamiltonian, V = a [X. s z ] d x , 
and to the order parameter, 5A = —A[X,s y ]. One 
can estimate the matrix elements of these operators as 
V nm ~ e(L z /l so ) and (SA) nm ~ A(L z /l so ), where e is 
the quasiparticle energy. Since we are primary interested 
in the effects of disorder, we will limit our discussion 
to the situation when the low-energy physics is disorder 
dominated, i.e. the random potential W > m&x{V, 5 A} 
masks the off-diagonal matrix elements of the determinis- 
tic Hamiltonian. Having this in mind we thus omit both 
V and 6 A throughout the paper. However, this assump- 
tion does not affect our results in qualitative ways. 



Assume that the chemical potential lies close to the 
bottom of the AT-th band, fi ~ fit • We refer to this 
band as the "topological band", since in the absence of 
interband scattering it defines whether the quantum wire 
is in the topological phase or not. The condition of the 
topologically non-trivial phase reads B 2 > B 2 = A 2 + 
(fi-fi z N) ) 2 . We also assume a hierarchy of energy scales 

lA N) - fif~ X) > E so > A ~ B, where E so = ma 2 /2 
is the spin-orbit energy. This condition implies that all 
other channels with the band index n < N are in the 
trivial phase, since for them B 2 -C A 2 + (fi — fti^) 2 - In 
the current experiments E so > A, i.e. our theory is on 
the border of applicability (see, however, the discussion 
in section IV ) . 



B. One-dimensional chiral fermions 

For strong spin-orbit coupling, E so 3> B, each band is 
characterized by two Fermi momenta, 



(11) 



as shown in Fig. [5j We aim to construct a low energy 
Hamiltonian describing the system at energy scales < 

E so . To this end we represent the fields ip^ as 

4 n \x) * R^(x)e ik > + L$(x)e-<*, (12) 
^\x) ~ L^{x)e- ik > + B$(x)e<*, (13) 

in terms of a superposition of right (R) and left (L) 
chiral fermions, and then linearize the Hamiltonian ( 10 ) 




FIG. 5. Chiral one-dimensional fermions in the spin-orbit 
quantum wire: a and b channels are shown by open and filled 
dots, resp. The magnetic field B affects the dispersion relation 
for the a-channel only if the chemical potential is close to /iz . 
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around the Fermi momenta. Modes with Fermi momenta 
of equal modulus define a conduction channel. We ob- 
serve that each channel belongs to one of the two subsets, 
a or b. In the a-channel i?-movers have spin up, and L- 
movers have spin down, while for the channels of type b 
spins are reversed. The a-channel of band index N plays 
a distinguished role in that it has zero Fermi momentum, 
k c ^ = 0. This particular mode defines what we call the 
"topological channel". For max{_B,A} <C E so it is the 
only channel strongly susceptible to the magnetic field, 
thus defining the topological phase of the whole wire. 



The effect of B (but not A) on other channels can be 
safely neglected. 

We next collect all R and L fields into two spinors, 

& = (-Raf) %i L a i, Lfrfi -Rat; Rbl, L a l, £&t)> (14) 

and W — (o'P h 5 , ) T where the band index (n) is left im- 
plicit. The spinor ^ acts in an 8 x TV-dimensional space, 
defined by the direct product of band (index n), chan- 
nel (a/b), chiral (R/L) and particle-hole spaces. In this 
representation the low energy Hamiltonian is given by 



where 



JV 



n = E 9 /* (n) (a0 



» 



-iA a RL <i 



(T ab S 



iA a RL 



v 

'<h {n) ) T 5 



~iv n d x - n (B/2) (t ab + of) S Nn 
(B/2) (t ab + af) 5 Nn iv n 8 x - ft 



■<& {m) (x)dx, 



RL 



(15) 



(16) 



Here. 



= [2m{ni N) - nT' + £ so )] 1/2 /m is the Fermi 
velocity of the n-th channel, the chemical potential is 
now defined relative to the energy /Uz , and W nm (x) is 
a random 8x8 matrix satisfying the hermiticity condi- 
tion (W nm y — W mn . In deriving this Hamiltonian we 
have neglected oscillatory terms with phases ix(±k a ±k{,). 
This is justified because the typical lengths involved in 
the problem (£ ~ A/v and Ib ~ B/a) exceed by far the 
Fermi length Xp ~ max(l/fc a , l/k ) determined by the 
large spin-orbit energy E so . We have also used that the 
order parameter A(x) can be chosen to be real. 
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C. Majorana representation 

The 8iV x 8N first-quantized matrix defining the 
Hamiltonian (151 satisfies to the p-h symmetry, 



(17) 



which is the defining condition for a class D Hamilto- 
nian (here, the transposition acts on kinetic term as 
d J = — <9 X .) However, all what follows will be more con- 



veniently formulated in an alternative representation, in 
which the symmetry assumes a different form: for each 
band n we define a set of eight Majorana fields 

= (i? a t + Rat)/V2, V R = (R at - R at )/V2i, 

g = (Rbi + Rbi)/V2, rtf = (R H - R H )/V2i, 

(and analogous relations for the L-movers, with spin re- 
versed), which we combine into the 8iV-spinor 

x= tf, g,£Z,v?, vtvtf- (is) 

The two spinors x and ^ are related by the unitary trans- 
formation 



X = U9, x T = *t/ f , 



(19) 



U = l ab 



1 RL 



V2 



1 1 

i —i 



We combine the Majorana fields of a and b channels as 
£ R = (£, R ,£, R ) T , and reorder the spinor components as, 

In this representation (which will be used throughout the 



rest of the paper) the Hamiltonian (151 takes the form 



N 1 f 

n,m—l 

H M= ( h^Snm + iW~ in a RL ® l ab 6 nm + iW~+ \ , . 

nm[ ) ~ ^ -in <? RL ® t ab 5 nm + w+- h^6 nm + w++ ) ' 

I 

where the deterministic part reads and the random matrices are constructed as 

ft<*> = -iv n a RL ®t ab d x - a RL ® , (22) W~ = ( "ST 1 , 

V W l W 2 J , 

A^=Aa?±(B/2)(l« b + *f)8 Nn , (23) f RR RL 

W- + = 'St W Al , (24) 
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in terms of real (symmetric) and imaginary (antisymmet- 
ric) parts of the random matrices (W) RL = w RL + i w RL 
etc. The remaining blocks of the IF-matrix are defined 
by 



..ab 



W++ = afW—af. 



(25) 

In the Majorana basis (20 1, the class D symmetry is ex- 
pressed through the antisymmetry H T = —H. 

We finally specify the statistics of disorder. We choose 
W(x) to be a (^-correlated and Gaussian distributed ran- 
dom matrix of size 8N x 8iV with a zero mean value 
(W(x)) and variance 



In 



6(x - x')(5u'8j>j> + SijfSjv), 



(w2(x)w2 J (a/)) = — S(x - x'){8 iV 5jj< - £y (26) 

Here, the composite indices (i,j etc.) label states in 
the direct product of band, channel and chiral spaces. 
The scattering matrices defined in this way break time- 
reversal and spin rotation symmetry (e.g., random spin- 
flip scattering caused by random spin-orbit terms is in- 
cluded in (26).) The strength of disorder set by the co- 
efficient j w translates into the "golden rule" scattering 
rate 



JV 



= 2 7 2£(1K) 



(27) 



of the normal conducting (i.e. 
pled) quantum wire. 



superconductor decou- 



IV. QUASICLASSICAL APPROACH 



The kinetic term in the low energy Hamiltonian (21) 



which was discussed in the previous section is linear in 
momentum, and this facilitates the formulation of quasi- 
classical equations of motion (aka Eilenberger equations) 
for the model at hancP. We here review the construc- 
tion of these equations in a manner closely following the 
spirit of Refs. [T^] and |2"U1 Af ter introducing the basics 
of the method (section IV A) wc construct the Eilen- 
berger Q-function in the limit of a single clean "topo- 
logical channel" (section |IV B[ ) and discuss the resulting 
density of states (section [IV C[ ). In section IV D we define 
the 1^2 topological invariant in terms of the Q-matrix. 
Section |IV E| outlines the general construction of the so- 
lution to the Eilenberger equation in the inhomogeneous 
disordered wire with boundary conditions. 



A. Eilenberger method 

We start by defining, 

x') = V^G*Li(x, x') o* L ^ 



where G R ^ A (x, x') = (x\(e±iQ—H)~ 1 \x l ) are the retarded 
and advanced Green's functions of the system. Under 
transposition (which in our current representation repre- 
sents the particle-hole symmetry) the function g behaves 
as 



i (x,x') = -af L [g^{x',x)] T a RL , 



(29) 



i.e. advanced and retarded Green's functions get inter- 
changed. It is not hard to derive the two mutually adjoint 
differential equations 



d x g e (x,x') +C e g e (x, x') = -i6(x - x'), 
d x > g e (x,x') - g e (x,x')£ e =iS(x- x'), 

describing the dynamical evolution of g. Here, 

C t = -i(u-V) 

where matrix Co has elements 



(30) 



(31) 



(32) 



and the operator V is related to the Hamiltonian matrix 
% as 



id x + v7 n 1/2 {af L H n 



-1/2 



(33) 



Due to the antisymmetry of % = — n 1 , the operator V 
obeys the particle-hole symmetry 



V = -a RL V T a? L . 



We next define the Eilenberger function as 

Q e (x) = lim \2i g e (x, x') — sgn(x — x')\ . 



(34) 



(35) 



where the subtraction of the sgn-function regularizes a 
discontinuity arising in g at x = x' due to the combina- 
tion of linear derivatives and (5-function inhomogeneity 
in Eqs. (30). Subtracting the two equations in (30), we 



then obtain the Eilenberger equation of motion 
d x Q e (x)+[L e ,Q e (x)] =0. 



(36) 



The Eilenberger function Q obeys the particle-hole sym- 
metry 



a? L Q T _ e (x)a 



RL 



-Qe(x), 



(37) 



and the normalization condition Ql{x) = 1, where 1 is 
the unit matrix (the latter condition can be checked by 
verifying that Eq. ( 36 ) preserve the normalization Q 2 = 



const.). The unit-value of the normalization constant is 
fixed by the jump height of the sgn-function in (35 1. We 



finally note that the operator V can be straightforwardly 
constructed from (33), however, for our present purposes, 



(28) we need not state its explicit form in generality. 



B. Eilenberger function in the clean limit 

As a warmup, we apply the quasiclassical approach 
to the limit (W — 0) of an infinite clean quantum wire 
subject to constant B, A. In this system all channels are 
decoupled, and we may concentrate on the 4x4 matrix 
Q e (B, n) describing the "topological" channel a with n = 
N. (The Eilenberger function of the other channels may 
be obtained by setting B = 0, rescaling the velocity and 
transforming A — > —A in the case of 6-type channels.) 

We start by introducing the 4x4 operator 



L. 



T RL 



-*A- 



r RL 



Iff 



RL 



e a 



(38) 



as the reduction of the general operator C € to a single 
channel. The solution Q e is then determined by the rela- 
tions [Q e , L e ] — and Q\ = 1. To solve these equations, 
we assume L e to be diagonalized as 



A = 



T(\< 
A+ 



JiL 



A_ 



(39) 



where the exact form of the (non-unitary) transformation 
matrix T will not be needed and 



A± = VAo ± A, 
A = B 2 + A 2 - M 2 -e 2 , 
A = 2 ^B 2 A 2 - (A 2 - e 2 



(40) 



are the eigenvalues. The defining equations for Q are 
then solved by matrices of the form 



Qe = TAT -1 , 



(41) 



where A = (±1, . . . , ±1) is a diagonal 4x4 matrix con- 
taining unit-modular entries in arbitrary configuration. 
The proper sign structure is determined by causality, i.e. 
the sign of the infinitesimal offset e — > e ± iO in the re- 
tarded/advanced Green function. That increment enters 
in the combination (e±iO)af' L , which means that the ap- 
propriate matrix structure of the retarded Green's func- 
tion (opposite for advanced) is given by 

A = <r* L . (42) 

A more explicit derivation of this structure is detailed in 
Appendix [A] 



C. Clean density of states 

The density of states in the bulk of the topological wire 
is given by v(e) — (2ttv)^ 1 Retr a RL Q € . The matrices T 
diagonalizing Q do not commute with af~ L , which means 
that a little extra work is required to evaluate the trace. 
We start from the representation 



p+ 



P~ 



(43) 



where P + and P are projectors on the space of L € 
eigenstates with eigenvalues ±A+ and ±A-, resp.: 



P+ =Tdiag(l 2 ,0)T- 1 , 
P~ =Tdiag(0,l 2 )T- 1 . 



(44) 



That this representation faithfully represents the matrix 
Q e is checked by application of (43) in the eigenbasis 



where all matrices assume a diagonal form. It remains 
to obtain a representation of P which does not make 
explicit reference to the diagonalizing matrices T. To this 
end, notice that (eigenrepresentation understood) L 2 — 
diag(A 2 f l 2 , A 2 ,l 2 ) = A 1 4 + AP+ - \P~. This equation 
can straightforwardly solved as 



P d 



Anl 



0-M. 



(45) 



Substituting this expression into ( 43 1 , we obtain a repre- 



sentation of Q which makes reference only to the operator 
(38), and its eigenvalues. Computing the trace, we ob- 
tain DoS profiles as shown in Fig. [6| Before discussing 
the structure of these results, a general remark may be in 
order: the DoS of a one dimensional quantum system is 
determined by an interplay of the kinetic energy operator 
(k f-> —id x ) and the "potential" (£). On the other hand, 
we computed the DoS from Q as determined by C, and 
this matrix seems to be oblivious to the kinetic energy. 
A closer look, however, shows that information on the 
band dispersion sneaks in via the nonlinear constraint 
Q 2 = t. Indeed, the conservation of the constraint, and 
the unit value of the normalization are consequences of 
the linearity of the derivative operator in ( 30 ) , which in 
this way co-determines the structure of Q. 

Inspection of Eq. (43) shows that the DoS contains 



singularities at the zeros A±(e) = 0, which arc located at 



e± = |S±v / A 1 



(46) 



Let us assume that B > A. Fig. shows the ensuing 
DoS profile, along with the underlying dispersion relation 
for a value of the chemical potential fi < fi c , where 



/'<-: 



y/B 2 - A 2 



(47) 



defines a critical value where the lower of the DoS sin- 
gularities, e_, touches zero and the band gap closes 
[Fig.[6]jb)]. At larger values [i > /i c , the gap reopens, (c), 
the DoS looks qualitatively similar to that of the /i < /i c 
regime, but the system is in a topologically distinct state 
(see the next section.) Finally, at values /x > fi*, where 



B 2 + B^B 2 + AA 2 /V2, 



(48) 



the lower band e_(fc) develops an extremum at finite val- 
ues of k which manifests in a third van-Hove singularity 
at the energy 



e - A^/l-B 2 /^ 2 , 
as shown in panel (d) of Fig. [6j 



(49) 
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(a) /j < (i c (b) /i = /U c (c) /i > /tic (d) /i > /i* 
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FIG. 6. Sketch of the density of states and the corresponding dispersion relations in the clean nanowire. The gap is closed at 
the transition point y, = fi c (a), and opened again (b,c). For fi > n* two minima develop (d). 



D. 



topological invariant 



The symmetry Eq. ( 37 1 implies that at zero energy 
the product af L Q e= o is an antisymmetric 4x4 matrix, 
which implies the existence of a Pfaffian. Due to the 
signature of A the determinant of Q is unity, the same 
is true for the determinant of the 4x4 matrix a^ L = 
(7% <S> la. Consequently, the Pfafhan of <r^ L Q e= o - 
which squares to the determinant of that matrix — may 
take one of two values, ±1. This motivates the definition 
of the topological index 



AT top = Pf Q e=0 ) = 



+1, 
-1, 



(50) 



distinguishing between the two phases. Computing iV top 
from (43), we find the index structure stated in (50). 
(Right at the critical point fj, — fi c , the matrix Q e= o 
becomes singular and the index cannot be defined.) 



E. Eilenberger equation with disorder 

In this section we discuss the formal solution of the 
Eilenberger equation in the presence of disorder. The so- 
lution is "formal" in the sense that the Eilenberger Green 
function will be a functional of a given realization of the 
disorder. To obtain practically useful information, one 
will want to average over different realizations, and this 
step of the computation needs to be done numerically, as 
discussed in the next section. 

To start with, consider the prototypical system geom- 
etry shown in Fig. [7] The terminals indicated at the left 
and right represent superconducting regions, assumed 
non-disordered for simplicity. (This is an inconsequen- 
tial assumption provided the rate of disorder scattering 
r^ 1 ru N"f1,{l/v n ) does not exceed the energy gaps (e_ 
or eo) in the terminals.) In these regions, the Eilenberger 
equation can be solved analytically, as discussed in the 
previous section. Describing the disorder present in the 
center region in terms of a generalized variant of the qua- 
siclassical evolution operator C, we will show how the left 
and right asymptotic configuration of the Green function 



get connected by a transfer matrix, M, functionally de- 
pending on the disorder configuration. The ensuing gen- 
eralized Green function will then be the starting point 
for our numerical analysis. 

To be more specific, we consider a quantum wire where 
the gap A(x) and/or chemical potential fx(x) vary in 
space in the region \x\ < L/2 and saturate to some con- 
stants A L / R and ^Iril a ^ x ^ —L/2 or x 3> L/2, respec- 
tively (Fig. [7]). These constants set asymptotic values of 
the Q-matrix, 



Q e {x -> -00) = Q. 



s(x 



-00) = Q+, (51) 



where Q± are constructed using the results of sec- 
tion |IVB| for the homogeneous profile of A, B and /Lt. 
The boundary Green's function Q_ and Q + may describe 
different or equivalent topological phases of the wire. 

We denote the Q-matrices obtained at the interface 
between the asymptotic superconducting regions, and the 
center disordered region, resp., as 



Qr = Q{xr), Ql = Q{xl), 



(52) 



where xr = —xl = L/2. These two configurations are 
related by a transfer matrix, 



Q R = M(x r ,x l )Q l M 1 (x r ,x l ). 



(53) 



For arbitrary positions x and x' the formal expression 
for the transfer matrix M(x, x') at given energy e follows 



Qi 



Q- 



XL 







—Of 







FIG. 7. Disordered spin-orbit wire connected to two ideal 
superconducting terminals which are described by the Eilen- 
berger functions Q+ and Q- . The Q-matrix at the boundaries 
between the scattering region and terminals is denoted by Qr 
and Ql- In the superconductors Q-matrix rapidly converges 
to either Q+ or Q_ on a scale of coherence length. 
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from the Eilenberger equation ( 36 1 



M e {x,x')=r x 



exp 




V{y)\ dy 



(54) 



with V x denoting the path-ordering operator. A relation 
similar to Eq. (53) connects the boundary matrices Qr/l 
with the Green's function in the far right/left region of 
the wire, 

Q± =M(x,x RjL )Q RIL M- 1 {x,x R i L ), (55) 

assuming that x — > ±oo. The transfer matrix satisfies 
certain symmetries. Along with the unitarity condition 
(cf. Eq. pi) 



o* L M]o* L =M7\ 



(56) 



it also satisfies to the p-h symmetry: with the use of 
Eq. ( 33 ) one obtains 



RL 



RL 



and this yields 



<j? L ML<J? L = M7\ 



(57) 



We now aim to represent the Q-matrix in the scattering 
region in terms of the transfer matrix M e and the asymp- 
totic Eilenberger functions Q±. We start by translating 



the transfer matrix relation ( 55 ) to a set of algebraic con- 



ditions relating the matrices Qr/l to Q±. To this end, 
notice that the action of the non-unitary (cf. Eq. (56)) 
transfer matrix on a generic matrix Qr will in general 
produce exponentially increasing and decreasing contri- 
butions. The former are inacceptable in that they lead 
to exponential divergencies in the quasiclassical Green 
function. As is detailed in Appendix [Bj the requirement 
of a non-divergent Green function leads to the algebraic 
conditions 



(1 + Q_)(l-Qi)=0, 
(1 + Qi)(l-Q_)=0, 



(58) 



while the right matrix Qr should obey the analogous 
relations 



(1-Q + )(1 + Q R ) = 0, 
(l-Q fl )(l + Q+) = 0. 



(59) 



We finally combined these equations with the transfer 
matrix relation ( 55 ) to obtain closed expressions for Ql/r 
in terms of the asymptotic configurations Q± and M. 
As a result of a straightforward calculation detailed in 
Appendix [B] we obtain 

Q * = 1+ Q ++ MQ.A/- 1 (1 ' Q+) ' (60) 

2 



1 + (1 



Q-+M-1Q+M' 



(61) 



where M = M(xr,xl)- 

These formulae define the starting point for an ef- 
ficient numerical computation of the disordered Eilen- 
berger function Q(x). To this end, one computes the 
transfer matrix M by numerical solution of correspond- 
ing system of linear first order differential equations. One 
next applies Eqs. (60) and (61) to obtain Qr/l- Finally, 
our main object of interest, Q(x), is obtained by appli- 
cation of M(x,Xrij j ) to either Qr or Ql- 

So far we have considered a quantum wire connected 
to two superconducting terminals. However, the gener- 
alization of the method to the system shown in Fig. [T] 
is straightforward. The key observation is that in the 
limit of vanishing barrier conductance jt«1 the chiral 
fermion fields satisfy 



R^{xl) = e l0 L bT (a; L ), Rbi(x L ) 



e^L ai (x L ) 



(62) 



where <f> is (energy dependent) reflection phase shift. 
(Here, we assumed the absence of barrier spin flip scatter- 
ing, inter-channel scattering or related complications.) In 
the limit of asymptotically high potential barrier <j) = tt. 
Relations ( 62 ) define boundary conditions for the Eilen- 



berger function Ql- As verified in Appendix D, these 
conditions assume the form of Eqs. (58), where, how- 
ever, the role of Q_ is taken by an "effective" matrix 
Q— = Q- (4>) describing the tunnel junction. The ex- 
plicit form of this matrix reads 



Q_ = (sin(<^K h + cos(0)a>f ) 



RL 



(63) 



This matrix also satisfies Q 2 _ — 1, and the rela- 
tions ( 60|61 1 stay intact. 

We close this section with an important statement: the 
Eilenberger functions Q+ and Q_ of the terminals define 
the "Majorana number'^ of the wire as 



M = Pf (<7 



RL , 



Pf(^Q-), 



(64) 



in terms of the product of two I12 topological invariants 
given by Eq. (50). It is shown in Appendix [Cj that for 
At = — 1 the system supports a Majorana fermion local- 
ized in the scattering region between two terminals. 



V. NUMERICS 

We next turn to the discussion of our numerical results 
obtained for the setup shown in Fig. [l] We have solved 
the quasiclassical equations according to the algorithm of 
section [IV E| and from there computed the local density 
of states (LDoS) v L (e) = (27T«) _1 Re tr(a^ L Q L ) at the 
left end of the wire close to the tunnel barrier. In the 
actual calculations wc shifted the energy into the com- 
plex plane, e — > e + iT. This shift accounts for the fact 
that in the "real" system states may escape to a con- 
tinuum of lead scattering states, which gives them the 
status of quantum resonances of finite lifetime ~ 
The corresponding decay rate^l is given by the standard 
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golden rule expression, T ~ grS where 5 ~ irv/2NL is 
the mean level spacing in the scattering region of size L. 
In principle, one might numerically compute the broad- 
ening by an extension of the numerical setup so as to 
include and extended normal metallic scattering region 
to the left of the tunnel barrier. This, however, would 
slow the performance of the numerics which is why we 
prefer to introduce the broadening "by hand" . At any 
rate, the Majorana peak, present for M. = — 1, will ac- 
quire a finite width ~ r. The DoS profiles obtained as a 
result of the procedure outlined above are presented and 
discussed in section HT1 above. 



Re(A±) > 0. This gives 

dp e ipx 



g(x,x';e) = T 



2lT p - i$ <g> a RL 



_ T -i = 



(A3) 



= - l -T{af L + sgn(x - x')) e^*^ T~\ 



Sending x' — > x and comparing to (35), we obtain the 
identifications (41) and (42). 



Appendix B: Boundary conditions 



VI. CONCLUSIONS 

In this paper we have adapted the Eilenberger quasi- 
classical approach to the specific conditions of a class D 
topological quantum wire. The most striking feature of 
the quasiclassical approach is that the Green functions 
of the system are described by ordinary, and hence eas- 
ily solvable differential equations. A numerical solution 
of these equations for given realizations of the disorder 
produces accurate information on the spectral proper- 
ties of reasonably realistic model systems hosting Ma- 
jorana boundary states. Our results confirm the predic- 
tions made in RefPfor a more idealized system, viz. that 
disorder generates spectral peaks which can be confused 
for genuine Majorana peaks. It thus seems that an unam- 
biguous detection of the Majorana might call for a mea- 
surement scheme beyond direct tunneling spectroscopy. 
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Appendix A: Derivation of Eq. (42 1 



The goal of this Appendix is to derive the structure 



( 42 ) by explicit calculation. Our starting point is the 



momentum representation of the equations of motion for 
the quasiclassical Green function g(x, x'; e) = g(x — x'; e), 



-p + -L e ) g(p,e) = 1. 
v 



(Al) 



Assuming a diagonal representation as in ( 39 ) , we obtain 

1 



9{P, e) = - T 



p - i- ® a RL 



(A2) 



The inspection of eigenvalues X± shows that for the re- 
tarded Green's function, i.e. at e — > e + iO, one has 



In this Appendix we derive the boundary condi- 
tions ( 58|59 l for the Q-matrix. The method we use is 
adapted from the circuit theory of RefsJMSS. p or sml _ 
plicity, we consider the 4x4 matrices describing individ- 
ual channels, generalization to many coupled channels is 
straightforward. 

Consider the quasiclassical function g(x, x') = 
g(x,x';e) to the right of the right terminal, x,x' > xr. 
By assumption, the generator L e 



X>X R 



is constant m 



space, and this means that the equations ( 30 ) defining 
g admit the solutions 

g(x,x R ) = e- L ^- x ^g(x R + 0,x R ), 
g(x R ,x) = g{x R ,x R + 0)e LAx - XR) , 

where we have again set v = 1 for notational simplicity. 
According to Eq. ( A3 1 , we have 



g{x R + 0,x R ) = --^{Qr 



g(x R ,x R + 0) 



■XQi 



1), 
1), 



(Bl) 



so that 



g(x,x R ) = --e- L ^-^(Q R + t), 

g(x R ,x) = -~(Q R + l)e L ^-*«l 

Defining g = T~ l gT and Q R — T~ 1 Q R T, we next trans- 
form to the representation ( 39 1 to arrive at 



g(x,x R ) 
g(x R ,x) 



- 2 (Qn 



t)e 



The key observation now is that the matrix functions 
g must remain finite as a; — > oo. Due do the posi- 
tivity of the matrix A, this is equivalent to the condi- 
tion (1 - a^ L )(Q R + 1) = (1 - af L ® l)(Q R + 1) = 
(1 — A)(Qr + 1) = 0. By the same token we have 
(Qr + 1)(1+ A) = 0. Finally, noting that TAT" 1 = Q+ 
represents the asymptotic form of the Q-matrix deep in 
the superconductor, we arrive at Eqs. (59). Eqs. (58) are 
shown in an analogous way. 
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We are now in position to derive Eqs. (601 and (611 



To this end, we take Eqs. (58) and multiply it by the 
transfer matrix M from the left and by M^ 1 from the 



right. Bearing in mind Eq. (55 1, one obtains 
1 - MQ_M~ l Q R + MQ-M~ X -Q+ = 



0. 



Adding this relation to Eqs. (59), we arrive at 
2 • 1 - MQ-Ar 1 Q R - Q+Q R + MQ-M~ l - Q H 



(B2) 



0, 
(B3) 



which in one more step yields the result (60). Eq. (61 1 is 
proven analogously. 



Appendix C: Majorana mode 

In this Appendix we discuss the "Majorana number" 
M. and show how the Majorana state emerges from the 
quasiclassical Eilenberger function. 

We consider the spectrum {Ej} of Andreev bound 
states which follows from the poles of the Q- 
function ( |60|61| ). If we denote 



23(c) = Q+(c) + M(e)Q_(e)M- 1 (c), (CI) 

then the energies Ej are solutions of the secular equa- 
tion defD(Ej) = 0. The Majorana state, if it exists, 
corresponds to Eq = 0. 

To proceed, we introduce matrices Q±(e) = Q±(e)<j]} h , 
and the secular matrix 23(e) = 23(e)<r RL . Since in the sub- 
gap interval of energies there is no distinction between the 
retarded and advanced Green's function, the unitarity re- 
lation {G R / A y = G A / R , the basic definitions Eqs. (2§ 
and ( [35] ) imply Q±(e) — —Q±(e). Further, particle-hole 
symmetry (37) yields e) = — Q±(e), Taking into 



account the class D symmetry of the transfer matrix, 



Eq. (57), the secular matrix takes a form 



23(e) = Q+(e) + Af(e)g_(e)M T (-e). 



(C2) 



It satisfies 23 T (e) = — 23(— e), and thereby guarantees that 
Andreev bound states appear in pairs ±Ej. 

One may ask now whether a zero energy solution of 
the secular equation exists or not. At e = 0, the particle 
hole symmetry puts tighter restrictions on the matrices 
Q± , M and 23 (we omit the energy argument for brevity) . 
One gets 



Q± 



? L M T a? L = M~ x 



Q±, 

M* 



M, 



23 T = 



-23, 23* = 23. 



(C3) 
(C4) 
(C5) 



We thus observe that both Q± and 23 are real antisym- 
metric matrices of size 87V x 8N, which enables us to 
rewrite the secular equation in terms of a Pfaffian 



Det23 



Pf(X3) = Pf(Q+ + MQ-M T ) 



0. (C6) 



Let us denote by fig the left null space of the matrix 23, 
i.e. any bra (0| £ Qq by definition satisfies (0|23 = 0. 
Since 23 is antisymmetric, the dimension of its null space 
is even, dim fig = 27V. At this stage we make use of a 
mathematical lemma proven in Appendix A of Ref. 1241 
the parity of the number Af can be expressed as 

= pf Q + pi [MQ- M T ] = DetMPf Q_Pf Q+. 

(C7) 



The particle-hole symmetry ( C4 ) implies that Det M = 
±1 at zero energy. Now, the transfer matrix M(x, x') is 
a continuous function of its arguments, with initial value 
M(x = x 1 ) = 1. We thus conclude that DetM = +1, 
so that the parity of Af is determined by the terminal 
configurations, 



(-!) Ar = PfQ_PfQ 4 



(C8) 



This parity is equal to the "Majorana number" M. intro- 
duced in section IIV El 

Let us now focus on the most interesting case M = 
1, corresponding, as we will see, to a single Majorana 
de^. In this case we have dim^Q = 2, and the null 



mo 



space is spanned by two linearly independent vectors. Let 
(<j>i I € fip b e the first basis vector. It is easy to check that 
(02 1 = (<j)i\Q + £ £Iq , and we may choose this state for 
the second basis vector in fig . Indeed, for any (</>i| e fig 
one has (4>i\T> — (4>i\T> — 0. Using the definition of 23, 



Eq. (CI), we deduce that 



(4> 1 \Q + M = -(fa\MQ-. 
This relation enables us to evaluate (02 1 23 as 
(0 2 |23 = (0 1 | (t + Q+MQ-M- 1 ) 

= (0i| + ((0i|Q+M) Q-M- 1 
= (0i| - (falMQ.Q-M- 1 = 0, 



(C9) 



(C10) 



and hence we proved that (02 1 € ^o- Since 23 is a 
real antisymmetric matrix, its right null space is ob- 

In other words, the two kcts 



tained as fi = (S1q) 



|0i,2) = ((0i,2)|)' satisfy the relation 23|0 12 ) = 0. 

Let us look at the Eilenberger function Q^(e) around 
its pole at e = 0. According to Eq. (60), it can be repre- 
sented by two equivalent equations, 



Q fl = l + 223- 1 (l-Q+), 
Q R = -1 + 2(1 + Q+)V-\ 

Multiplication by af 1 " yields 



Q^3 RL = ^3 RL + 223- 1 (l-Q+), 
of h QRaf h = -a 3 RL + 2(1 - Ql)V-\ 



(Cll) 
(C12) 

(C13) 
(C14) 



At e —¥ the inverse operator from the secular matrix 
has the pole structure, 23 _1 ~ 2?./e, where the matrix 1Z 
is its residue at zero energy. 

At this stage it is advantageous to introduce a new bra 



(X±| = (0l|±(0 2 | 



l|(i±Q+), 



(C15) 
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and ket 



\x±) = &±Q+)\<h), (C16) 



basis in the null spaces JIq and flff, resp. The bra basis 
(x±| i s n °t orthogonal and we can denote by \r]±) the ket 
basis, which is dual to it, (Xa\v<r') = &aa'- Using these 
definitions, we may formulate resolutions of unity, 

1 = l»?<r)0c<r| = \Xa){v<y\, 

where a = ± is summed over. 

Now let us look at the singular part of the matrix Qr, 

-K(l-Q+), 



(i - Ql)n. 



(C17) 



We note that the matrix P_ = §(1 — Q+) acts as the 
projector in the space S7q , since 



(X+|P_=0, and ( X -\P- = (x-\- 
Similar relations hold in the ket space fi^: 

pT| x+ ) =0 , and P*\ X -) = \X-). 
Equivalently, we may write 



(C18) 
(C19) 



P_ 



\v-Kx-\, 



pi 



Ix-X^-I- 



Using these properties, it follows a^ L Q R n6 a^ L = 
-1ZP- — -P^TZ. Multiplying these relations by P_ from 
the right, and using the projector property P_ = P_, we 
obtain the relation 



RL rising ML 



a z Q 



-pIkp- 



\ X -)(V-\K\ri-}(x-\, 



or 



4"P__ 



*? L \x-)(x-W 



RL 

z ) 



where we defined 1Z = (r]—\TZ\rj—). Taking into ac- 
count our definition of the quasiclassical Green's func- 
tion, Eqs. (28 1 and (35), we finally obtain that the 



singular part of the propagator around zero energy at 
x = x' = xr takes the form 

G(x R ,x R ;e) ~ v-^ L \ X -)(x^ L ^ 1/2 , 

(C20) 

where v is the diagonal velocity matrix. 

This expression can be compared with the spectral de- 
composition of the Green's function, 



G(x, x , e) 



\Mx))(Mx')\ 



dp 



E 



dp 



\^ Ej (x))(^ Ej (x')\ 



E 4 



\^e(x))(^e(x')\ 
e-E 



where \i/je{%)) are normalized eigenfunctions of the BdG 
Hamiltonian, E g is the gap in the spectrum of the wire, 
the sum is going over the set of subgap Andreev levels 
(we have particle-hole symmetry P_j = —Ej) and the 
first term is present if the system contains a Majorana 
state. One thus concludes that 



1 



-1/2 



RL 



X-) 



\^o{xr)) 



(C21) 



is the amplitude of the Majorana particle at the point 
Xr. The amplitude of the Majorana state at any other 
point x can be obtained from \i}) (xr}} by applying the 
transfer matrix M(x, xr). 

To conclude, we have shown that if the Majorana 
number of our system is topologically non-trivial, i.e. 
M. = — 1, then the Green's function has the pole at zero 
energy. Up to a prefactor, the residue at zero energy is 
then the projector on the one dimensional linear subspace 
spanned by the Majorana particle. 



Appendix D: Matrix Q_ of a tunnel junction 

In this appendix we show that the boundary condi- 
tions for the Eilenberger matrix Q of a spin-orbit wire 
terminated at the left end by a (infinitely high) tunnel 
barrier are equivalent to algebraic relations ( 58 ) with the 
effective matrix Q_ given by Eq. 



Let us work in the original particle-hole basis where 
the spinor ^ takes the form specified by Eq. (14). We 
start by rewriting the left boundary conditions ( 62 ) in a 
matrix form. If one defines the reflection matrix 











(Dl) 



in the particle space, then f* is the reflection matrix for 
holes. Consequently, the spinor ^>(xl) satisfies the con- 
dition 



l-RW(x L )=0, R = 



f 
f* 



(D2) 



The full reflection matrix R is unitary obeying the 
particle-hole symmetry, o-P h P T crP h = R. Hence the 
boundary condition for the bar spinor ^ = (<7 P ! h, I') T takes 
the similar form 



(D3) 



*(x L )[l-R) =0 



According to the definition (28), the Green's function 



g e (x,x') inherits the boundary condition right to the 
tunnel junction from those of the direct product of two 
spinors, *&(x) ® {^>(x') a^ L ) . We thus conclude that 

(t-R)g(x L ,x;e) = 0, 
g(x,x L ;e)(t + R) =0. 
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When deriving the second relation, we have taken into 
account that af L r af L = — f. The above conditions are 
valid for any x > xt,. Sending now x — > Xl + and using 
the relations (Bl I written for the matrix Ql, one obtains 



(1-R)(t-Q L ) = 0, 
(l + Q L )(l-R) = 0. 

We see that these boundary conditions are equivalent to 
the algebraic conditions (58 1 if one identifies Q_ = — R. 



The normalization R 2 = 1 guaranties that Q- belongs to 
the manifold of the quasiclassical Eilenberger functions. 
Transforming matrix Q- into the Majorana representa- 
tion (20) we finally obtain the result (63). 



The scattering phase 4> is the parameter of the matrix 
Q_ which may depend on the band index n and charac- 
terizes the tunnel junction. In our numerical simulations 
we have used <p = it for all bands, which corresponds to 
the infinitely high barrier as compared to the energies of 
Andreev bound states. 
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